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We present semiclassical approximations to Green’s functions of multidimensional systems, extending 
Gutzwiller’s work to the classically forbidden region. Based on steepest-descent integrals over these functions, 
we derive an instanton method for computing the rate of nonadiabatic reactions, such as electron transfer, 
in the weak-coupling limit, where Fermi’s golden-rule can be employed. This generalizes Marcus theory to 
systems for which the environment free-energy curves are not harmonic and where nuclear tunnelling plays 
a role. The derivation avoids using the ImF method or short-time approximations to real-time correlation 
functions. A clear physical interpretation of the nuclear tunnelling processes involved in an electron-transfer 
reaction is thus provided. In the following paper, ^ we discuss numerical evaluation of the formulae. 


I. INTRODUCTION 

Electron transfer is a key step in many important 
molecular processes, including redox reactions in elec¬ 
trochemistry and charge separation in photosynthesis 
and solar cells.The electron resides initially on a 
donor molecule and is transferred to an acceptor, ac¬ 
companied by a reorganization of the polar environ¬ 
ment. This reaction can be characterized as a transi¬ 
tion between the donor and acceptor electronic states 
with potential-energy surfaces describing the reactant 
and product environments. As such, electron transfer is 
the simplest example of a nonadiabatic reaction involving 
transitions between different electronic states, requiring 
a theoretical treatment beyond the Born-Oppenheimer 
approximation.We are thus interested in studying a 
multidimensional curve-crossing problem, which as it in¬ 
volves discrete states, is inherently quantum mechanical. 
Although we pay particular attention to the electron- 
transfer problem, the nonadiabatic formalism is also rel¬ 
evant in other areas of science, such as defect tunnelling 
in solids.® Most formulae derived in this paper can be 
transferred directly to these fields. 

In this paper, we consider the golden-rule limit which 
occurs when the transfer of the electron, rather than the 
rearrangement of the environment, is the bottleneck to 
the reaction.^ This is quite commonly the case in prob¬ 
lems of interest, especially if the donor and acceptor are 
spaced far apart. The standard approach for treating 
electron-transfer problems in this limit is Marcus theory,® 
which is based on Fermi’s golden rule with the additional 
approximation that the nuclei are treated classically. The 
free-energy curves of the environment in the reactant and 
product states are also assumed to be harmonic with the 
same frequency. Although it is quite common for the 
environment fluctuations to be harmonic,there is 
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no reason for the frequencies to be exactly equal in the 
two cases, unless the reaction is symmetric such that the 
products are equivalent to the reactants. In some cases 
the reorganization energies for the forward and backward 
reaction can differ by a factor of two,^® and it is neces¬ 
sary to use a more general rate expression which allows 
for this asymmetry but retains the classical and harmonic 
approximations.^®’^"^ 

Marcus theory also neglects nuclear quantum effects, 
such as tunnelling and delocalization, which have been 
found to be significant in electron-transfer problems even 
at room temperature.®® It is however possible to com¬ 
pute the quantum golden-rule rate exactly for the spin- 
boson model,®®’®' which treats all environment modes 
as linearly-coupled harmonic oscillators. Extensions of 
this to treat non-linear couplings is also possible.®® The 
multilayer multiconfigurational time-dependent Hartree 
(MCTDH) method®® is in principle able to compute the 
exact rate for such systems,®® but in practical applica¬ 
tions to large systems is usually limited to specific forms 
of the Hamiltonian such as system-bath models. 

Thus for many problems in the golden-rule limit, an 
accurate calculation of the reaction rate will require an 
efficient method which includes nuclear quantum effects 
and treats the potential-energy surfaces in a general 
way without making global harmonic approximations. 
In this article, we present a new derivation of such a 
method for computing the rate approximately using a 
time-independent picture for the nuclear degrees of free¬ 
dom. 

Our derivation is based on an exact expression for the 
golden-rule rate in terms of the Green’s functions de¬ 
scribing the nuclear quantum dynamics of the reactant 
and product systems at a given energy. By extending 
Gutzwiller’s work,®®^^® we present the semiclassical limit 
of these Green’s functions in the classically forbidden re¬ 
gion, which may also be useful for other applications. 
These functions are defined in terms of imaginary-time 
classical trajectories which join together to form a peri¬ 
odic orbit, also known as an instanton. The definition 
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of the rate is found by performing a number of steepest- 
descent integrations. 

Other instanton approaches are well known from adia¬ 
batic rate^^^^® and tunnelling splitting^^’^^ calculations, 
where in both cases the Born-Oppenheimer approxima¬ 
tion is first applied to obtain a single-surface Hamilto¬ 
nian. The Imi^ method^^ can be used to derive the in¬ 
stanton approximation to the adiabatic escape rate from 
metastable states in the high- or low-temperature limit^® 
although its application to finite temperatures^® or nona- 
diabatic reactions®® is less well understood. In this ap¬ 
proach, the rate is expressed approximately as the imag¬ 
inary part of a free-energy, which is obtained by analytic 
continuation of a divergent integral.®^ 

Although instanton theory dates back many years, a 
numerical method for its efficient application to com¬ 
plex multidimensional systems has only been developed 
more recently,^®’®® and it is thus experiencing a revival 
of interest.®®^®® The present work can be considered to 
belong to the same class of methods, which treat the col¬ 
lective motion of the nuclear modes with a multidimen¬ 
sional instanton. We note however that our approach is 
very different from the semiclassical instanton approach 
presented in Ref. 40. This was derived to recover Mar¬ 
cus theory in both the normal and inverted regimes using 
an instanton to describe the dynamics of the transferred 
electron, while treating the environment classically. 

It has become common to define thermal reaction rates 
using time correlation functions.^® However, our deriva¬ 
tion in energy-space gives us access not only to the ther¬ 
mal rate but also the energy-dependent reaction proba¬ 
bility. This can be combined with any energy distribu¬ 
tion (including microcanonical) to give a nonequilibrium 
instanton rate, which may be helpful for understanding 
scattering calculations of large molecules or gas-phase 
unimolecular reactions. The modified correlation func¬ 
tion presented in Ref. 42, which is used for computing 
nonadiabatic rates avoiding oscillatory functions, also re¬ 
quires that the flux operator is biased towards energy- 
conserving electron-transfers. If we are to combine this 
method with an instanton approach, it will be natural to 
use a formulation in energy-space rather than in time to 
access the necessary variables. 

An outline of the paper is as follows. We express the 
formula for the golden-rule rate in terms of Green’s func¬ 
tions in Sec. H and give a semiclassical approximation 
to these functions in Sec. HI. Using this approximation, 
the golden-rule instanton rate is derived in Sec. IV in 
two forms which are evaluated analytically for a system 
with linear potentials, the spin-boson model, and in the 
classical limit for general potentials in Sec. V. Sec. VI 
concludes the article. In Paper H,® which follows this ar¬ 
ticle, we discuss how the instanton formula can be applied 
numerically to complex systems and relate it to Wolynes’ 
quantum instanton approach.^® 


II. QUANTUM GOLDEN-RULE RATE 


We consider a general multidimensional system with 
two electronic states, each with a nuclear Hamiltonian of 
the form 


1=1 

where n € {0,1} is the electronic-state index and 
X = (xi, ... ,Xf) are the Cartesian coordinates of / nu¬ 
clear degrees of freedom. These nuclei move on the 
potential-energy surface 14,(x) with conjugate momenta 
Pj = Without loss of generality, the nuclear de¬ 

grees of freedom have been mass-weighted such that each 
has the same mass, m. The electronic states |n) are cou¬ 
pled by A(x) to give the total Hamiltonian in the diabatic 
representation,® 

i? = i?o|0)(0|+i7i|l)(l|-f A(x)(|0)(l| + |l)(0|). (2) 

The thermal rate, k, of transfer of reactants, defined 
by the projection operator |0)(0|, to products, defined by 
|1)(1|, is given exactly by 

kZo= P{E)e-^^dE, (3) 


where /? = l/k^T is the reciprocal temperature and Zq 
the reactant partition function Tr [e~^^ |0)(0| ] • PiE) is 
the reaction probability with energy E, defined as^® 


P(U) = i(27r?i)®Tr FS{E - H)FS{E - H) 


(4) 


where the flux from reactants to products is 

F=^A(x)(|0)(l|-|l)(0|). (5) 

Note that, although we shall not make use of it in this 
paper, it would be simple to replace the canonical ensem¬ 
ble in Eq. (3) with any distribution of energy to obtain 
microcanonical or non-thermal reaction rates. 

In this paper, we consider only systems for which the 
electronic coupling, A(x), is very weak such that second- 
order perturbation theory, known as the golden-rule ap¬ 
proach, can be employed. In this limit the formulae re¬ 
duce to 


P{E) = Tr A{x}6{E - Ho)A{x)S{E - Hi) 


Zn = Tr 




( 6 ) 

(7) 


These equations, combined with Eq. (3), give the golden- 
rule rate in the form also derived from a limit of the 
non-oscillatory flux correlation function,^® and are exact 
to order A®. 
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Expanding the trace in Eq. (6) in the basis of reactant 
and product internal (e.g. vibrational) states |p,) and |i/) 
would give the standard golden-rule formula, 

^^0 = T E I (^1 Im ) I' - E .) > (8) 

fiiy 

where the sum over states should be replaced by the inte¬ 
gral /f dEfj^dEi, for continuous systems. Considering this 
golden-rule formula, we notice that the reaction occurs 
between internal states of equal energy. For the major¬ 
ity of systems of interest (in the Marcus normal regime), 
at low enough temperatures the dominant contribution 
comes from low energy states which overlap only in the 
classically forbidden region. In this case therefore the re¬ 
action includes a nuclear tunnelling process and requires 
the treatment of such quantum effects to be described 
adequately. 

However, unless the internal states are known, which 
would be the case for instance if the potential-energy sur¬ 
faces are assumed to be harmonic, this formulation can¬ 
not be applied straightforwardly. In the present work, 
we seek a semiclassical approximation which allows us to 
evaluate the rate in the golden-rule limit for complex sys¬ 
tems by avoiding the computation of the wave functions 
explicitly. 

We introduce the Green’s functions (fixed-energy prop¬ 
agators) describing nuclear dynamics on each of the elec¬ 
tronic states independently, defined equivalently by the 
following two expressions 

G„(E;) = lim (E + ir;-i/„)-i (9) 

ri^0+ 

= lim QiiE+ivWh 

^^ 0 + n Jo 

The imaginary part is related to the density of states by 
ImG„(i?) = —Tr6{E — Hn). The reaction probability can 
thus be written as 

P{E)=d yy'A(x')Im(x'|Go(E)|x") 

X A(x")Im(x"|Gi(£;)|x')dx'dx". (II) 

A similar formula was derived by Miller et al.^^ for adia¬ 
batic reactions in which it was also noticed that only the 
imaginary part of the Green’s function need be known in 
order to obtain the reaction rate. 

An exact evaluation of the Green’s functions will be 
impossible for complex systems as this is equivalent to a 
complete solution of the Schrodinger equation. Nor does 
the microcanonical density operator treated in a path- 
integral representation"^®’^® lead to a practical computa¬ 
tional technique. We therefore require a simpler semi¬ 
classical description of the imaginary part of the Green’s 
functions in the classically forbidden region. The deriva¬ 
tion of this is given in the following section. 


III. SEMICLASSICAL GREEN’S FUNCTIONS 


Like the wave function which solves the Schrodinger 
equation for the Hamiltonian in Eq. (1), the Green’s func¬ 
tion defined by Eq. (10) contains all information required 
to study the nuclear dynamics. It would therefore be a 
very useful object to obtain and apply to a wide range of 
problems, although it is in general as difficult to compute 
exactly as the wave function itself. However, it can easily 
be defined using Feynman’s path-integral propagator, 
from which one can take a semiclassical approximation 
replacing the path integral as a sum over classical trajec¬ 
tories. 

In many previous applications, semiclassical Green’s 
functions were employed to describe quantization in 
bound states.This required locating periodic orbits 
in the classically allowed region, for which there may be 
numerous possibilities, many of which are unstable, espe¬ 
cially in large complex systems exhibiting chaotic dynam¬ 
ics. However, in this work, the most important quantum 
effect is that of tunnelling and we neglect quantization 
in the reactant and product potential wells. This means 
that we are interested in evaluating the Green’s func¬ 
tion only in the classically forbidden region. As we shall 
show, this requires us to locate a single trajectory, and 
is therefore expected to lead to computationally feasible 
methods even in complex systems. 

Gutzwiller^®^^^ has derived a semiclassical approxima¬ 
tion to the Green’s function, (x'|G„(E)|x"), in the classi¬ 
cally allowed region, i.e. where H„(x') < E and 14(x") < 
E. The derivation starts from the van-Vleck propagator, 
a semiclassical approximation to (x'|e“‘'^’*^* 


iL„(x',x",t'-n= E 

cl. traj. 


ASr,/h 

(27ri7i)-f 


G„ 


dx'dx" 


where the full action is 


( 12 ) 

(13) 


5„ = ^„(x',x",t'-0 


r' 

1 

dx{t) 

L 

2TO 

dt 


K(x(t)) 


dt. 


(14) 


Here, x(t) is a classical trajectory which travels from x" 
at time t" to x' at time and we sum over all such 
trajectories. To avoid clutter, the transpose is implied on 
the second partial derivative variable of Hessian matrices 
throughout, e.g. ■ As the Hamiltonian is time- 

independent, we can without loss of generality set t” = 0. 
We note that — gives the energy conserved along the 
trajectory. 

Using van-Vleck’s propagator in Eq. (10), Gutzwiller 
employed a stationary-phase evaluation of the Laplace 
integral. The stationary-phase points solve 

E5„(x',x",4) + E = 0, (15) 

utn 
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which picks out times, tn, defining trajectories with the 
required energy. This gives a semiclassical approximation 
for the Green’s function as a sum over classical trajecto¬ 
ries with energy 


G„(x',x",£;)= ^ 

cl. traj. 


AWn I'h—IV'KI2 /-jC'i 

(2m?i)(/+i)/2 ^ ^ ’ 


Now the same classical trajectories can are formally de¬ 
fined as continuous paths starting at x" and ending 
at x' and giving a stationary value of the abbreviated 
action, 


Mq)=y! 

Wn = Wn{x',x",E)= pr^dq, (17) 

Jx(q)=y." 

which is a line integral along the trajectory x{q) repre¬ 
sented by the generalized coordinate q. The magnitude 
of the momentum at any point along the trajectory is 

Pn = \/2m[E - t4(x)]. (18) 


The prefactor is found from the determinant of an 
(/ -|- l)-dimensional square matrix. 


= (- 1)^+1 


d'^Wr, d'^Wr, 
ax'Ox" d^'dE 

dEdx" dE^ 


(19) 


Following Gutzwiller^^, we take the absolute value of this 
determinant and introduce a phase term in Eq. (16) de¬ 
termined by the Maslov-Morse index, v, which is given 
by the number of conjugate points along the trajectory. 
For our purposes, conjugate points occur where the tra¬ 
jectory encounters a turning point, where Ki(x) = E, and 
comes instantaneously to rest.^^ We call this a bounce. 

We now derive an equivalent semiclassical Green’s 
function formalism for the classically forbidden region, 
where 14,(x') > E and Vn{x") > E. Some work in this 
direction has already been completed to describe tun¬ 
nelling using imaginary-time trajectories, using similar 
approaches to that outlined in this paper. Some general¬ 
ization is however required as results in the literature con¬ 
sidered only pathways which pass through the forbidden 
region but with end points outside,or treated only 
certain one-dimensional systems.We also introduce 
a transition-state theory approximation below which al¬ 
lows us to present a simple expression for the Green’s 
function between classically forbidden end points in the 
general multidimensional case. 

Using the semiclassical van-Vleck propagator in 
Eq. (10) gives stationary-phase points which again solve 
Eq. (15). As before, solutions correspond to classical tra¬ 
jectories travelling from x" to x' in a given time. However, 
classical dynamics in the forbidden region is only possible 
in imaginary time, as it is well-known that this is equiva¬ 
lent to real-time dynamics on an upside-down potential- 
energy surface.Therefore, whereas in the classically al¬ 
lowed region these stationary points lay on the real-time 



FIG. 1. Representation of the deformed contour for the 
Laplace integral in the classically forbidden region. The two 
marked stationary-phase points correspond to the direct and 
bounce trajectories. 


axis, for the forbidden region they are on the imaginary- 
time axis. The action for these trajectories is purely 
imaginary. 

Let us consider the nature of these imaginary-time tra¬ 
jectories. In the /-dimensional coordinate space, there 
will be a (/ — l)-dimensional “turning” surface on which 
Vn{x) = E. This surface is a generalization of the turn¬ 
ing point, which separates the allowed from the forbid¬ 
den region, and imaginary-time trajectories can bounce 
off it. Possible solutions therefore correspond to trajec¬ 
tories which travel directly from one end point to the 
other or those which bounce of the surface in-between. 
Because the dynamics are time reversible, a bouncing 
trajectory will always travel along the same pathway 
before and after the bounce, which always occurs nor¬ 
mal to the turning surface. The direct trajectory corre¬ 
sponds to td = —Ud) and the bounce to % = —irb, where 
R) > Td > 0. The stationary point of the former is a max¬ 
imum along the imaginary axis but a minimum along the 
real axis, whereas the latter is vice versa. Note that the 
reverse of these trajectories also exist with Imt > 0 but 
occur in the wrong part of the complex plane to be of 
interest for the contour integration which we shall per¬ 
form. 

For the simplest cases of tunnelling into an unbound 
potential, such as the one-dimensional linear Vn{x) = 
KnX, only these two classical trajectories can exist at a 
given energy. However, a simple one-dimensional bound 
potential, such as the harmonic oscillator, also sup¬ 
ports trajectories which pass through the turning sur¬ 
face into the classically allowed region where it bounces 
any odd number of times, before returning to the forbid¬ 
den region.®^ These correspond to complex times with 
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Ret > 0 and contribute to the real part of the action 
and hence to the phase of the Green’s func¬ 
tion. However, we are ultimately interested in condensed- 
phase systems where the potentials are multidimensional 
and give rise to chaotic dynamics.This causes the 
phases to mostly cancel each other out when we take 
further integrals over the Green’s functions. We shall 
therefore assume that we can ignore the stationary points 
with Re t > 0 and compute a steepest-descent integration 
along the contour 7^ -I-7b -I-700 shown in Fig. 1 such that 
it passes through the two stationary points. A small pos¬ 
itive value of 77 causes the function to decay with Re t and 
ensures that the integral along 700 is 0 . 

In this way, we are restricting the information needed 
about the system to the classically forbidden region. Our 
approximation is a form of quantum TST as we are ef¬ 
fectively ignoring coherence effects at long times. The 
coherences are responsible for the discrete energy levels 
of a finite bound system,and the zero-point energy of 
an infinite system.In many condensed-phase systems, 
these effects do not play an important role and can be 
safely ignored. Only at very low energies will we find a 
serious error due to the discrete nature of states in the 
wells which cannot be described using this approach. 

More complex potential forms may support many 
bounce or direct trajectories as well as those which 
bounce more than once. In such cases, one can either 
sum over all possibilities or turn to a path-integral sam¬ 
pling scheme which avoids the steepest-descent integra¬ 
tion in position variables. In this work, we derive only 
a semiclassical steepest-descent rate theory and shall as¬ 
sume the simplest case that only one trajectory of each 
type exists, but still without specifying the exact form 
of the surface. We shall show in Paper iR how Wolynes’ 
quantum instanton method,which can be used for 
systems with rough potentials, is related to our golden- 
rule instanton approach. 

We can then follow a very similar derivation to that of 
Gutzwiller’s in the allowed region performing steepest- 
descent integrations about the points td/b- One must 
take particular note that our integration along 7 b passes 
over only half of the Gaussian peak centred at tb- This 
gives the following semiclassical approximation to the 
Green’s function: 


g„(x',x",f;) = rd(x',x", A) + ir^(x',x", A), (20) 


where the functions are defined either by the di¬ 

rect or bounce imaginary-time trajectory as appropriate. 
These trajectories are continuous paths starting at x" and 
ending at x' which give a stationary value of the imagi¬ 


nary abbreviated action. 


= w„(x',x",f;) = 



Pn = \/2m[Vni'X.) - E] 


r^/"(x',x",if) 


D„ 


2TT^/\Dn\ -W,^/n-iu7r/2 

(27r^)(/+i)/2 


(- 1 ) 


/-tl 




Ox'ax" 


dx'QE 


dEdx" 


dE2 


( 21 ) 

( 22 ) 

(23) 

(24) 


We use the notation of a bar over all variables related 
to imaginary-time propagation. Note that the imaginary 
action W is so called because it is defined as —iW. In 
this region, like pn = —ipn, it is always real and positive. 
Unlike the Green’s functions in the classically allowed 
region, which were oscillatory, these decay exponentially 
with W. 

The direct trajectory, with i/ = 0, therefore contributes 
to the real part, and the bounce, with zx = 1, to the imag¬ 
inary part of the semiclassical Green’s function. Even if 
trajectories exist with more than one bounce, they can 
normally be ignored as their large actions will ensure that 
they do not dominate either the real or imaginary parts. 

The semiclassical formulae suffer in the same way as 
the Wentzel-Kramers-Brillouin (WKB) approximation^^ 
from poles at the turning points of a trajectory. That 
is, at a turning point where Un(x) = E, the prefactor 
Dn goes to infinity, as can be seen by transforming 
to a coordinate basis parallel and perpendicular to the 
trajectory.^^ However, we have found them to be a good 
approximation to the exact Green’s functions when x' and 
x" are far from the turning points, which as we shall show 
is all that is required for a steepest-descent evaluation of 
the golden-rule rate constant. 

We have shown that the Green’s function is not just a 
simple analytic continuation of Gutzwiller’s formula be¬ 
cause a factor of a half appears in the imaginary part. 
This is in agreement with previous work®^’^^ which ex¬ 
plain the one-dimensional WKB connection formulae. 
Interestingly the factor also appears in the Im F approach 
where it is argued that only half of the imaginary barrier 
partition function is required.The derivation given 
here based on the contour integral used to compute the 
Laplace transform of the van-Vleck propagator is more 
direct and requires no analytic continuation. 

For illustration, we shall give a simple example of the 
semiclassical Green’s function for the one-dimensional 
linear potential-energy surface Vn{x) = KnX. The semi¬ 
classical approximation in the classically allowed region, 
Knx\ Knx" < E, is given by 


I 


Gnix',x'\E) = -- 




Pn{x’)pn{x”) . 




( 25 ) 
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with 


= 


1^1’ = 


Pn{x'' f - Pn{x')^ 


3mKn 

Pn{x”)^ + Pn{x')^ 

|3tok„| 


(26a) 

(26b) 


or in the forbidden region, K„a;', Knx" > E, 

i 


Gn{x',x",E) = 




‘h V Pn{x')pn{x") 


+ 2^' 


(27) 


with defined as Eqs. (26) with bars added over 

the variables and momenta defined as in Eq. (18) and 
Eq. (22). This is the same semiclassical result as found 
in Ref. 54. 

The wave functions for this linear potential are 
known to be ^/>„(a;;i?) = s/a Ai{z{x; E)), where a = 
and z{x;E) = (2m/K^?i^)^/^(K„a; — E). 
The exact Green’s function for this linear potential can 
be found using the procedure outlined in Ref. 58 from 
the Wronskian, w. It is 



z{x;E) 

FIG. 2. The diagonal elements of the Green’s function for 
a linear potential and plotted with real and imaginary parts 
in red and black against a dimensionless function of position 
and energy. The exact form, Eq. (28), is shown with solid lines 
compared to the semiclassical approximations, Eqs. (25) and 
(27), with dashed lines. 


(x'|G„(£;)|x") = ^Ai(z>)Ai 




w = m 


— 1 + i/ 


(28) 

(29) 


The usual procedure for performing steepest-descent 
integrals assumes that pre-exponential terms are approx¬ 
imately constant over the range in which the exponential 
dominates. That is for functions A(q) and B{q) of a 
d-dimensional variable q. 


where z< and are the lesser and greater of z{x'\E) 
and z{x"-, E). 

The semiclassical approximation is compared to the 
exact Green’s functions in Fig. 2, showing good agree¬ 
ment everywhere except near the classical turning point 
z = 0. In fact the imaginary part of the semiclassi¬ 
cal Green’s function approaches the exact result asymp¬ 
totically as can be proved using Im (a;'|G„(A)|x") = 
—TTipnix'; E)il)n{x"\ E) and an asymptotic approximation 
for the Airy function,®® accurate for | 2 ;| —> oo given by. 


Ai(z) 


exp[- 


2VEz1/4 

sin[§(-z)3/2-K 


z>Q 
z < 0. 


(30) 


Note again that the factor of half appearing in Eq. (20) 
was necessary for this equivalence. 


IV. GOLDEN-RULE INSTANTON FORMULATION 

In this section, we shall derive an approximate formula 
for the golden-rule rate using two semiclassical Green’s 
functions in the classically forbidden region. We shall 
perform this derivation using a steepest-descent integral 
first over the positions and then over energy in the exact 
expression for the rate. This picks out two imaginary- 
time classical trajectories, which when joined together 
are known as the instanton. 



dq = (27r;i)‘^/®B(qt) 


a® A ”2 

q=qt 


(31) 


where ql is defined such that A(ql) is a minimum. This is 
also known as the semiclassical approximation because, 
as long as B{q^) ^ 0, it gives the term with lowest order 
of h correctly. The error in the approximation is always 
an order of h higher and thus becomes exact if —?► 0. 

This steepest-descent approach requires that only one 
minimum of the function A(q) exists, but is easily ex¬ 
tended to treat multiple well-separated minima by sum¬ 
ming over the contributions. However, like other in¬ 
stanton approaches,the resulting formulation will 
not be able to treat the rough potential surfaces found 
for reactions in liquids where many local minima exist 
whose steepest-descent integrands would overlap. Such 
problems are better treated using path-integral Monte 
Garlo or molecular dynamics approaches,which 
we discuss in the following paper. ^ The instanton ap¬ 
proach derived here is only directly applicable to solid® 
or gas-phase systems as well as system-bath models of 
condensed-phase electron transfer.However, it is the 
derivation and physical interpretation of a rate formula 
which is the focus of the present work, and so we will as¬ 
sume for now that the potentials are sufficiently smooth. 

As already pointed out above, for the calculation of 
electron-transfer rates it is the classically forbidden re- 



















7 


gion which dominates the result. We note that the 
Green’s function will be much simpler to work with nu¬ 
merically in this region, where except for phases orig¬ 
inating from the conjugate points, there is a real ex¬ 
ponent, and it is thus not an oscillatory function. By 
a similar principle other imaginary-time path-integral 
calculations®^ are numerically tractable, including, for in¬ 
stance, instanton calculations of adiabatic reaction rates 
and ring-polymer transition-state theory (RPTST).^®’®® 

In Sec. IV B we transform our rate expression from the 
language of the Hamilton-Jacobi formalism, where the 
trajectories are defined by their energy, to the language 
of Lagrangian dynamics, where they are defined by the 
elapsed time. 

As our formula for the rate will be derived from 
a steepest-descent integration over the position coordi¬ 
nates, a consistent calculation of the reactant partition 
function is 


1 

^0 « TT -^, (32) 

2sinh \l3fiujj 

where iOj are the normal-mode frequencies at the min¬ 
imum of Vo(x)- This form assumes that there are no 
translation or rotational degrees of freedom but can be 
easily modified to treat gas-phase problems by project¬ 
ing out such modes in the usual manner. In this case the 
steepest-descent approximation is equivalent to a local 
harmonic approximation for the region near the bottom 
of the reactant well. However, this is only a minor ap¬ 
proximation, and in what follows, we do not assume that 
this harmonic approximation holds for the whole reactant 
potential. 


A. Hamilton-Jacobi formalism 


We insert our semiclassical approximation for the 
Green’s functions into Eq. (11). This requires only the 
imaginary parts which are, according to our approximate 
formula, Eq. (20), given by trajectories with a conjugate 
point, i.e. those which bounce. Note also that the ener¬ 
gies of the Green’s functions on each electronic state are 
equal. We thus obtain the semiclassical reaction proba¬ 
bility 


Psc{E) 



A(x')A(x") 

(27r?i)/-i 


^-w/h (jx'dx". 


(33) 


where the combined action is IT = Wq + Wi and the 
trajectories considered are those which bounce once be¬ 
tween the end points x' and x" as shown in Fig. 3. We 
have assumed that there is only one such trajectory on 
each surface. If more than one exists, one should either 
take only the dominant trajectory with the smallest value 
of W or sum over all possibilities. Note that the integral 



FIG. 3. Schematic showing the two imaginary-time bounce 
trajectories with energy E near the crossing point in a one¬ 
dimensional, two-state system. The blue trajectory is on |0) 
and the red on |1). The steepest-descent integration of posi¬ 
tions will be taken about the crossing point x' = x” = x^ at 
which Vo{x^) = Vi(a:*) = V*. 


in Eq. (33) shall be computed using the method of steep¬ 
est descent which avoids including the spurious poles of 
Dn. 

The steepest-descent integration is taken about the 
point x' = x" = xl-, at which 

dW 

— = p' -p; = 0 (34a) 

dW 

^ = + = (34b) 

where (or its equivalent with double primes) is the 
imaginary momentum of the trajectory on surface n at x'; 
it has magnitude Pn(x') and direction pointing along the 
trajectory. As classical mechanics is time-reversible, the 
direction chosen is immaterial, but should be consistent. 

The fact that the energy E is equal for both trajecto¬ 
ries by construction implies that lies on the crossing 
seam. This has physical significance showing that al¬ 
though the bouncing instanton trajectory is delocalized, 
the hop between electronic states occurs dominantly on 
the crossing seam. For a one-dimensional system this 
point is uniquely defined, but in general, the seam is 
an (/ — l)-dimensional surface and the hopping point x^ 
varies for trajectories of different energy. Because the 
imaginary momenta on each surface cancel at the hop¬ 
ping point, the trajectories must join smoothly into each 
other to form a periodic orbit. This, along with the con¬ 
straint that both trajectories must reach a turning point 
defines the hopping point x^ for energy E. Note how¬ 
ever that this does not imply that the momentum at this 
point is necessarily normal to the Vb(x) = Vi(x) surface. 

In the following, all terms shall be evaluated at this 
hopping point x' = x" = x^. This includes the electronic 
coupling for which we therefore need only one value, de¬ 
noted A = A(xt). The semiclassical reaction probability 









is thus defined as 


1 

(35) 

This is the form given for a system in the so-called 
normal regime where the reactant and product minima 
lie on opposite sides of the crossing seam. When they lie 
on the same side, known as the Marcus inverted regime, 
it is clear that a different Ansatz would be required to 
define the instanton because Eqs. (34) no longer dehne 
the stationary point. 

The result in Eq. (35) could be used directly to com¬ 
pute microcanonical reaction rates.The approximation 
is only valid for energies lower than the activation energy 
and will also deviate strongly from the exact result 
at very low energies in condensed-phase systems. This 
is due to the transition-state theory assumption which 
ignores the true density of states in the reactant well 
and neglects its zero-point energy. We concentrate how¬ 
ever on computing the thermal rate which is found by 
integration over all energies weighted by the Boltzmann 
distribution. It would also be possible to extend the the¬ 
ory to use other energy distributions and describe certain 
nonequilibrium effects. 

Inserting the approximation for the reaction probabil¬ 
ity into Eq. (3) and performing a steepest descent integral 
over E gives the semiclassical result that we seek: 


Stt r ~^—^— 

Psc{E) = —A'^^/DoDi 


d^W 

d^W 


d^W 

dx'dx" 

d^W 


dx"dx' dx"dx" 


kscZo = V2TTh ^ v DqDi 


d^W 

dx'ax’ 

d^W 


d^W 

dx'dx" 

d^W 


dx"dx' dx"dx" 
1 


d^W 

dJ^ 


^-W/h-/3E 


(36) 


where the value is given at energy E which solves -I- 
fdh = 0 for a given temperature. The full derivative im¬ 
plies that x' and x" are not held hxed but are moved to 
the appropriate value of x^, which is the stationary point 
of W for each given energy. Of course in one-dimensional 
systems where the hopping point is always the same, 
there is no difference here between the full and partial 
derivative with respect to E. 

Because t„ = — gives the imaginary time of the 
trajectory, the energy is chosen by the steepest-descent 
integration is that which ensures that the total imaginary 
time taken by the orbit is Tq-I-ti = ph. We have therefore 
obtained a formula similar to the semiclassical descrip¬ 
tion of the quantum Boltzmann distribution,®® which also 
leads to the usual derivations of adiabatic instantons in 
terms of imaginary-time periodic orbits.The golden- 
rule instanton is thus a periodic orbit of constant energy 
and total imaginary time PU. It follows a continuous path 
on the Vb(x) surface from x^ to a turning point before re¬ 
tracing its steps back to x^. Here it hops to the Vi(x) 
surface, without needing to modify its momentum as the 



FIG. 4. Schematic showing three instanton paths (blue and 
red) in a two-dimensional system at different temperatures. 
The shortest path corresponds to the highest temperature 
which tunnels at high energy, whereas at low temperature, 
the tunnelling pathway is longer and at low energy. Con¬ 
tours are shown for both potential-energy surfaces such that 
one contour level is chosen to be the energy E of each of 
the instantons. In each case, this contour level surrounds the 
classically allowed region where the imaginary-time trajectory 
cannot enter and the instanton bounces normal to this sur¬ 
face. The hop occurs at x* without affecting the momentum 
of the orbit and is represented by the change in colour of the 
pathway. These points vary with temperature but are always 
located on the crossing seam Vb(x) = Vi(x), shown by the 
dotted line. 


potentials are equal, and performs a similar bounce be¬ 
fore returning. The periodic orbit must retrace itself after 
the bounce because it must approach and depart from a 
turning point, at which it comes to rest, along the poten¬ 
tial gradient. In Fig. 4, we show how the periodic orbit 
and the hopping point x^ are affected by the temperature 
parameter p. Higher temperatures lead to shorter in¬ 
stantons and thus a less delocalized and more “classical” 
reaction where tunnelling is less pronounced. Lower tem¬ 
peratures allow more freedom for the instanton pathway 
to become curved, an effect known as corner-cutting.^^’®® 
Equation (36) is only one of many possible ways to 
formulate the golden-rule instanton rate and we consider 
now transforming it to an alternative form. Applying the 
chain rule at the stationary point defined by Eqs. (34) 
gives 


d^IT d'^W 52 IT dx' 52 IT dx" 
dE2 “ dEdx' dA dEdx" dE 

and 


(37) 



d 

dE 


dW\ 

( d^W 

ax' 1 — 

dx'dE 

dW 1 

1 d^W 

dx"/ 

\dx"dE 


/ d'^W d^W 


dx'dx' 

d^W 


dx'dx" 

d^W 


\dx"dx' dx"dx" 



( 38 ) 
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Solving these linear equations leads to the following relationship between full and partial derivatives: 


d^W 

dJ^ 






( a'^w 

a^w ^ 

\ ^ / a^w 

1 a^w 

a'^w \ 

ax'ax' 

1 a^w 

ax'„ax" 


V dEax' 

aEax" ) 

a^w 



\ ax" ax' 

9x"9x" y 

' \ax"aE 


(39) 


such that the semiclassical rate can be written equivalently as 


kscZo = V2Trh ^\JDqDi 


d^w 

a>cax' 

d^W 


d^W 

dx'dx" 

d^W 


9 x'' 9 x' 

a^w 

dEdx' 


d^W 

dx'd_E 

a^w 


ax''ax" ax''dE 
a^w a^w 
dEax" aE^ 


1 

2 


^-w/n-pE^ 


(40) 


which is of course what would have been found by performing the steepest-descent integral over positions and energy 
in one step. 


Equation (36), and equivalently (40), which give a 
semiclassical approximation to the golden-rule rate is the 
main result of this work. If the instanton trajectory can 
be located and derivatives of its action evaluated, the 
expression can be applied directly to complex systems. 
We describe numerical methods for doing this in Paper 
II based on optimizing discrete pathways. 

We consider our derivation to be simpler and more 
rigorous than previous golden-rule instanton approaches 
based on an extension of the ImE approach.This pro¬ 
cedure was applied by Cao and Voth^^ to describe nona- 
diabatic transitions, although there appears to be no 
physical argument to explain the use of the Im F premise 
in this case. Here no saddle point was found in the spa¬ 
tial degrees of freedom but instead the imaginary part 
came from an analytic continuation of the divergent in¬ 
tegral along an imaginary time coordinate. We shall how¬ 
ever find that the formulae they obtained is very similar 
to that derived in this work. This shows more clearly 
that these instanton approaches are a form of transition- 
state theory (TST), i.e. real-time dynamical effects are 
ignored. This would be harder to understand from the 
ImE approach where time does not appear. 

Our derivation was performed using the Hamilton- 
Jacobi formulation to define the classical trajectories, 
which was the natural choice for treating trajectories 
which were required to have equal energies. However, 
there also exists an alternative formulation of classical 
dynamics using the Lagrangian picture, i.e. with a given 
imaginary time rather than energy. In the next section 
we perform a Legendre transformation to find an equiv¬ 
alent definition of the semiclassical rate formula. 


action®^ 

5 „ = 5 „(x',x",t„) = 


r 

hm 

9x(r) 

Jo 


dr 


+ 14(x(r)) 


= IT„(x',x",£;) + £;t„, 


dr 

■(41) 

(42) 


where the trajectory, x(r), travels through the classically 
forbidden region from x(0) = x" to x(t„) = x' with en¬ 
ergy E = ^1^. Because tq + ti = pti, the exponent of 
Eq. (40) is simply S = So + Si, which is the total action 
of the periodic orbit. Because the energies of the two 
trajectories are equal by construction, our values for t„ 
must also obey or = 0 where tq = Ph — r, 

Ti = r and thus . The hopping points are 

as before x' = x" = xb Note that such trajectories with 
equal end points are forced to bounce as long as t„ is not 
zero, i.e. 0 < t < pti. 

Equation (42) is a Legendre transformation and ex¬ 
presses the relationship between the two equivalent dy¬ 
namical formalisms. Eollowing Ref. 34 and taking deriva¬ 
tives, we obtain the relations 


d'^Wn 

dE^ 

d^Wn 

dx'dE 

dx'dx" 


(d^s^y^ 

V drl ) 

d^Sr, (d^sy-^ 

dx'dTr, 5r2 J 

d^Sn _ d^Sn f^y' d^Sn 
dx'dx" dx'dTn V J 9t„9x"’ 


(43a) 

(43b) 

(43c) 


B. Lagrangian formalism 

The Lagrangian formulation of classical mechanics 
is based on Hamilton’s principle function, Eq. (14). 
The imaginary time-version of this gives the Euclidean 


and the equivalents with any exchange of x' and x". 

Therefore, if the derivatives of Sn are known, it is 
a simple matter to identify all derivatives of IT„. Us¬ 
ing Eqs. (43) and standard properties of determinants, 
Eq. (Al), some simple but laborious algebra gives 
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and 


= (- 1)^+1 

' d^S, 

r - 


d^W„ 

dE^ 

-1 

a 


52IT„ 

d^Wn 

fd^Wn\ 

' d^Wn 

dx'dx" 

dx'dE 

{ dE^ J 

dEdx" 


dx'dx" 


where 


E = 


02^0 92^1 


a^s 

dx'dx" 

_ d^S 

dx"dx' dx"ax" 


d^S 

dx'dx' 

a'^s 


a^s 


a^s 


dTdx' drax" 


E 

n,n'—0 


a^s 

dx'dT 

a^s 


ax''dr 

a^s 

dr^ 


a^w 

a" IV 

a'^w 

dx'dx' 

a^w 

dx'dx" 

a^w 

dx'dE 

a^w 

Ox" ax' 

a^w 

ax''ax" 

a^w 

dx'idE 

a^w 

dEdx' 

aEdx" 

dE^ 

a^s„ 

dx^dx' 

d^Sn 

9^S„ ^ 

dx'dx" 

a^s^ 

1 ^nn 

ax" ax' 

dx" dx" / 

1 


d^So d^Si 


5r2 


i9r^2 


a s„ 

dx'dTr, 

a^s„ 

. dx" drn 


and Sn'n' = 2(5nn^ — 1 is a permutation symbol taking values ±1. 


/d^So d^Si 




9 S^, 

aT„,ax' 


9 

aT„,dx" 


(44) 

(45) 

(46) 

(47) 

(48) 

(49) 


We are now able to re-express Eq. (40) as 


kscZo = V27r?i 




s/n 


(50) 


or following a similar transformation to that of Eq. (39), 
with W and E replaced by S and r, as 


kscZo 




dr^ J 


C = 


a^s 


dx'dx' 

a^s 


ax"ax' 


a^s 

dx'dx" 

a^s 


ax"ax" 


' e-^/^ 


(51) 

(52) 


This is very similar to the golden-rule instanton formula¬ 
tion of Cao and Voth^^ which was derived from an exten¬ 
sion to Im F theory. The only difference is that an extra 
approximation to the prefactor was made, valid only for 
the spin-boson model (see Sec. VB), that the determi¬ 
nants cancel with the partition function to give 

_ 1 

fcc&v = (53) 


where again r is given as the value which solves ^ = 0. 
However, had this approximation not been taken, the two 
approaches would have given exactly equivalent results. 
It is interesting to see that the current approach gives the 
same result as an Im F formulation which was applied to 


nonadiabatic problems without rigorous derivation. The 
work presented in this paper, which gives an equivalent 
result, thus provides an explanation of why the golden- 
rule instanton formulation of Cao and Voth recovered 
the semiclassical result for the spin-boson model,and 
would also apply to more general systems if the extra 
approximation had not been made. In fact, this shows 
that the Im F approach works remarkably well in the 
golden-rule limit. However, the new derivation based on 
semiclassical Green’s functions is simpler as it does not 
involve analytic continuation of divergent integrals and 
starts from a rigorous energy-space picture of the reac¬ 
tion. 


Generalizations of this nonadiabatic instanton ap¬ 
proach have been proposed which interpolate the 
electron-transfer rate between the weak- and strong¬ 
coupling limits.®®’®^ However, they were also based on 
extensions of the Im F approach, which does not appear 
to lead to a unique formulation. In one case, the in¬ 
stanton was projected onto pure diabatic states, which 
causes errors in the adiabatic limit,®® whereas in the 
other, all electronic configurations were allowed giving 
a mean-field approach which would fail to describe the 
high-temperature golden-rule rate.®^ when the instanton 
shrinks. 
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V. ANALYTIC SOLUTION IN SPECIAL CASES 

It is possible to find a closed-form expression of the 
golden-rule instanton rate in a few special cases. In par¬ 
ticular, we treat a system with linear potentials, the 
spin-boson model and the classical limit for a general 
condensed-phase electron transfer. 


A. One-dimensional linear potentials 

The simplest description of a nonadiabatic curve cross¬ 
ing problem is that of two linear potentials, 

V„{x)=V^+KnX (54) 

with Ato > 0 > Ki and a constant coupling, A. This is 
similar to the famous Landau-Zener model®® but with¬ 
out the constraint that the particle travels at a constant 
speed along the x coordinate. In order to define a rate, we 
assume that the potential only has this linear form near 
the crossing point and flattens out at extreme values of 
X without affecting the transmission probability for rele¬ 
vant values of energy. The reactant partition function Zq 
could then be defined using the translational invariance 
of incoming scattering states. 

Using the wave functions given in Sec. Ill, and noting 
that the integrals over the Airy functions can be per¬ 
formed analytically using Ref. 69, gives the exact trans¬ 
mission probability 


P{E) = 


—4m^ 


1/3 


X Ai 


h Ko{ko — Ki)^Ki 

/2to(ko — ^ 

V Kqk\ 


E 


■ (55) 


The golden-rule instanton approach, Eq. (35), gives for 
E < 0, 


Psc{E) 


ttA^ I 2m 

7i|Ko — KiI V Vi — E 

- E)V^\ko - Ki\ 
X exp --———^-• 


(56) 


Using the asymptotic limit, Eq. (30), shows that this 
tends to the exact result as E ^ —oo. It diverges how¬ 
ever ai E = Vi and is not applicable for higher energies. 
Eor the thermal rate, defined by Eq. (50), we obtain 


ksG 


fcci exp 


24m{Ko — Ki)^ ’ 


where the classical rate is^^ 


(57) 



C(A - Vt) 



FIG. 5. Transmission coefficient and thermal rate constants 
for the one-dimensional linear model. Solid lines are exact 
results and dashed from instanton approximation. Energies 
and temperatures are weighted by ^ = {2mj to give 
dimensionless units. 


We plot these results in Eig. 5 for the symmetrical case 
K = kq = —Ki- As expected the semiclassical results 
deviate from the exact transmission probability when 
the energy nears the diabatic crossing Vi. Interest¬ 
ingly however, there appears to be a cancellation of er¬ 
rors as the thermal rate gives the true value perfectly. 
This is confirmed formally using the integral identity 
J^^Ai^{—aE)e~^^dE = j for a,/3 > 0. 

The same result is of course obtained using the golden- 
rule flux correlation function^^’^^ and the exact path- 
integral propagator. 

The figure shows that tunnelling changes the rate by 
orders of magnitude at low temperature but that the 
classical result is correct in the high-temperature limit. 
This simple example could be used to give a rough es¬ 
timate of the importance of nuclear tunnelling in more 
complex systems if a one-dimensional reaction coordinate 
is known. However this is not generally the case and the 
full-dimensional instanton pathway will normally need to 
be found for accurate predictions. 


B. Spin-boson model 


kciZo = 


1 27rm A^ 
j3h'^ h\Ko - Ki| 




The standard model for electron-transfer reactions in 
(58) the condensed phase is the spin-boson model.The 
potentials of reactants and products are given by sets of 
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shifted harmonic oscillators, 


/ 


Vb(x) = \muj]{xj 

3=1 

(59a) 

f 

^i(x) = Y - G 

1=1 

(59b) 


where r can be chosen to simplify the numerical integral 
over t as much as possible but has no effect on the result. 

A common approximation to the quantum golden-rule 
rate takes a stationary-phase integral about the point 
t = —IT such that 0'(—ir) = 0. This gives^^ 


and the electronic coupling, A, is taken to be a constant. 
The frequencies and couplings of the modes are deter¬ 
mined according to a given spectral density, which we 

write the discretized form, J(w) = | ~ 

where Cj = 

The exact quantum golden-rule rate for this system 
can be calculated numerically using^' 


«QM 


A 2 /-oo —ir 


^ J — cio —i- 


= —let -\— J 


1 — cos UJt 


uj'^ [tanh \l3fiuj 


-I- i sin wt 


(60) 

dw, 

(61) 


fcsp 


Ti^ 





(62) 


Similar formulations exist for generalizations of the spin- 
boson model to include anharmonicities,^^ although prac¬ 
tical expressions are limited to system-bath problems of 
a particular form where the bath can be represented by 
an effective harmonic form.^^ 

We now proceed to compute the golden-rule instanton 
rate for this system. The Euclidean actions for classical 
trajectories in this harmonic system are^^’^^ 


/ 

5'o(x',x",ro) = ^ 

3 = 1 

5i(x',x",ri) 

3 = 1 


rriLUj 

2 sinh ujjTq 
mujj 

2 sinh uJjTi 


[((2;^ + ^3^ + + ^ 3 )^) coshojjTo - 2{xt -f + ^j)] 

[{{x'j - + ix” - coshwjTi - 2{x'j - ^j){x'- - - eri. 


(63a) 

(63b) 


Solving 11^ = = 0, we find that 



-0 


sinh ^ujj{l3h — 2r) 
sinh Uj 


(64) 


where Uj = f3hujj/2. Therefore the periodic orbit which 
hops at this point has action 


The second derivatives are 


d^Sn 

dxtdx'f. 

d^Sn 

dxtdx'l 


= 8 . 


'jk 


mojj 


tanha;,T„ 


= -8 


'jk 


mu3j 


sinhwjTn 


(66a) 

(66b) 


and their equivalents with all single and double primes 
exchanged. This gives 


/ 

S'(r) = —er + ^ 2mu3j^^ 
3=1 


1 — coshwjT 

-5- - -hsmhw.-T , 

tanh Uj J 

(65) 


which is simply <()(—ir) for an /-dimensional system. We 
must in general choose r numerically to solve = 0 
to ensure that the energy of the two Green’s functions 
are equal. In the unbiased case where e = 0, symmetry 
considerations give r = ^f3h. However, we continue with 
the proof in the general case. 


c„ = n 


sinhd 


'3 ‘n 


j 

_ CoCi 

72 ’ 


2-f 


tanh ^ujjTq 


tanh ^ujjTi 


tanh ^uijTi \ 

tanh ^uijTo J 


(67) 

( 68 ) 

(69) 


where we have used Eq. (A3) as well as a number of 
trigonometric relations. For this case of the spin-boson 
model, the determinants in Eq. (51) cancel with the par¬ 
tition function, Eq. (32), and we may proceed with the 
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formula of Cao and Voth, Eq. (53), without approxima¬ 
tion. It is then seen that our instanton approach gives 
the same result as the stationary-phase approximation, 
Eq. (62), for the spin-boson model. Note that for this 
harmonic system, the fluctuations about the instanton 
pathway are treated exactly and so the same result would 
also be obtained by the quantum instanton method of 
Wolynes.^^ 

It was shown for example in Ref. 15 that this 
stationary-phase approximation introduces an error of 
about 20% for a model of aqueous ferrous-ferric elec¬ 
tron transfer. Such an error is quite acceptable as often 
only the order of magnitude of the rate is required. The 
stationary-phase approximation, and thus the golden- 
rule instanton, is not equivalent to the “semiclassical” 
method described in Ref. 73. In the latter case, an extra 
approximation was made that Hq and Hi commute and 
thus gives inferior results at low temperatures.^® 

The stationary-phase approximation to the golden- 
rule rate in the spin-boson model becomes exact in the 
classical (high-temperature) limit and recovers Marcus 
theory.®’^^ It is interesting to note that adiabatic instan- 
tons tend to an asymptotic low-temperature limit but 
break down at a certain cross-over temperature and a dif¬ 
ferent form is needed in the high temperature regime. 
This does not occur in the instanton approximation for 
the golden-rule rate, not just for the spin-boson model, 
but for any system as we shall show in the next section, 
which deals with the high-temperature limit more gener¬ 
ally. 


C. Classical limit 

We now consider a limit which has not been previously 
well studied, which is the classical limit for the golden- 
rule rate of a generic anharmonic system. We recently 
derived a general formula for the classical golden-rule 
rate^^ and we will show here that the semiclassical in¬ 
stanton method tends to a steepest-descent version of it. 

To achieve the classical, high-temperature limit, we 
let /? —>■ 0 which forces the periodic orbit to become an 
infinitesimally short line. As before we require that the 
hop occurs at x' = x" = x^ where the potentials are equal, 
and the trajectories must bounce along the direction of 
the gradient and be continuous at the hop. It is therefore 
possible to conclude that at the hopping point, x^, the 
directions of the two gradients are exactly opposite and it 
is thus the minimum on the crossing seam, Vo(x) = Vi(x). 
We are therefore allowed to describe the potentials using 
the following series expansion about this point: 

I4(x) « 0 -b g(l^(x - x*) -b ^m(x - x*)^Q„(x - X*). 

(70) 

We choose an orthogonal coordinate system such that the 
degree of freedom j = 1 is normal to the crossing seam 
at the hopping point, x^. All other modes are at their 
minimum positions here such that g„ = (k„, 0,..., 0). 


Note that this local expansion is by no means equivalent 
to the global harmonic approximation employed by the 
spin-boson model and Marcus theory. Here for instance, 
the frequencies at the crossing point are allowed to dif¬ 
fer on each electronic surface as well as in the reactant 
well. As the instanton method is derived by perform¬ 
ing steepest-descent integrations, including higher order 
terms in the expansion would not change the following 
results. We thus expect to obtain the correct activation 
energy even when the reorganization energy of the reac¬ 
tants and products are unequal. 

In the short-time limit, 

5„(x',x",t„) = ^\ x ' - x"p -b 14 (Kx' -bx")) T„, (71) 

and we find that in order for the total action to be sta¬ 
tionary with respect to x^, 


To = /3h 


Kl 

Kl - Kq ’ 


Ti = 


KQ 

Kq- Kl' 


The fluctuation terms tend to 


Cn = 



(72) 


(73) 


and due to symmetry between x' and x" we can simplify 
the determinant by rotating the axes to (x" -bx') /^/2 and 
(x" — x!)/i/2. This is equivalent to using Eq. (A2), and 
gives 
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(74) 


(75) 

(76) 


where f2 = and the minor |Q|ii is formed by 

removing the first row and column of the matrix and 
taking the determinant; it is defined to be equal to 1 in 
the case that / = 1. Therefore the golden-rule instanton 
method approaches the following form in the classical 
limit: 


1 27rm 


A2 Zt 


fol.hTST ^ 


- 0 Vt 




(77) 

(78) 


For a one-dimensional system, = 1 and the for¬ 
mula is obviously equal to the classical golden-rule TST 
rate.i’^^’^i Equation (77) is proportional to this one¬ 
dimensional formula evaluated along a particular reac¬ 
tion coordinate, which is defined as being normal to the 
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crossing seam. This rate is simply multiplied by classical 
fluctuations (vibrational partition functions in the har¬ 
monic approximation) along perpendicular coordinates, 
as expressed by the minor of the matrix. We name this 
harmonic classical golden-rule TST and note that it re¬ 
covers Marcus theory when applied to the spin-boson 
model. It therefore makes the same assumption that 
the reaction coordinate can be separated as is used in 
Born-Oppenheimer classical TST,^® and is therefore the 
golden-rule equivalent of this. A simple extension would 
give the equivalent to Eyring TST^® by replacing the clas¬ 
sical partition functions in by their quantum versions. 
However such a method is not as accurate as the instan- 
ton approach as it still treats motion along the reaction 
coordinate with classical dynamics. 

It is a major success of the golden-rule instanton the¬ 
ory presented here that we recover the classical limit for 
high temperatures. This is however not unexpected as 
we have started from exact golden-rule expression and 
performed a semiclassical approximation thus preserv¬ 
ing the term of lowest order in U, which of course gives 
the correct steepest-descent classical result. Semiclassi¬ 
cal instanton methods developed to study single-surface 
reactions^^’^® break down at a particular crossover tem¬ 
perature where the instanton collapses to a singularity 
and require a different formula for the high-temperature 
limit. This occurs because the potential at the barrier 
tends to a parabola with finite curvature which can¬ 
not support short periodic orbits above a certain cross¬ 
over temperature.However, the golden-rule instan¬ 
ton never collapses completely because the potentials be¬ 
come approximately linear near the activation energy 
forming a cusp, and at least a small amount of tunnelling 
occurs at all temperatures. 

We also note that the extra approximation taken by 
Cao and Voth,®® gives a method, Eq. (53), which assumes 
that the normal modes at the minimum of the reactant 
potential are equal to those at the transition state on each 
surface. This therefore only gives the correct classical 
limit for the case of the spin-boson model. 


VI. CONCLUSIONS 

In this paper we have derived a semiclassical instanton 
method for computing the rate of an electron-transfer 
reaction in the nonadiabatic limit. Our derivation starts 
from an exact expression for the golden-rule rate written 
in terms of Green’s functions, which are themselves ap¬ 
proximated by a semiclassical limit based on classical tra¬ 
jectories. The remaining integrals are evaluated within 
the steepest-descent method. This procedure defines two 
imaginary-time trajectories which contribute to the rate, 
one on each of the reactant and product potential-energy 
surfaces, both of which have the same energy and must 
encounter a bounce. The hop between the surfaces oc¬ 
curs at a point where the potentials are equal and the 
trajectories join together smoothly to make a periodic 


instanton orbit with imaginary time PTi. 

In this way, we have derived four equivalent formulae 
to define the golden-rule instanton rate: Eqs. (36), (40), 
(50) and (51), which use either a Hamilton-Jacobi or La- 
grangian formulation and either full or partial derivatives 
with respect to energy or imaginary time. In Paper H, 
which follows this article, we shall show how these for¬ 
mulae can be evaluated numerically for a complex mul¬ 
tidimensional system. 

In our opinion this derivation makes clear the assump¬ 
tions and approximations being made and provides phys¬ 
ical insight into the electron-transfer process, showing 
that the dominant contribution comes from an instan¬ 
ton of constant energy which hops on the diabatic cross¬ 
ing seam. Previous golden-rule instanton formulae were 
based on the Im F approach®® and it is interesting to see 
how closely the methods relate to each other, considering 
that their derivations are apparently so different. 

The golden-rule instanton approach gives the exact 
thermal rate for a system with linear potentials and, 
for the case of the spin-boson model, reproduces the 
well-known stationary-phase evaluation of the quantum 
golden-rule formula, Eq. (62). Note that although the 
steepest-descent approach uses only harmonic fluctua¬ 
tions about the dominant instanton pathway, the impor¬ 
tant integral along the instanton path is evaluated ex¬ 
actly. The remaining approximation is equivalent only 
to a local rather than a global harmonic approximation. 
In other words the harmonic frequencies are allowed to 
vary along the instanton pathway and the method should 
also give good results for anharmonic systems, including 
problems where the reorganization energies of products 
and reactants are different. 

In the high-temperature classical limit, it reduces in 
general to a steepest-descent form of classical golden- 
rule TST^^ and thus to Marcus theory for the spin-boson 
model. It can therefore be seen as an extension of clas¬ 
sical golden-rule TST to the quantum regime. How¬ 
ever, unlike the adiabatic instantons there is some tun¬ 
nelling at all temperatures and no crossover from the 
deep to shallow-tunnelling regimes. This shows that nu¬ 
clear quantum effects are always apparent when comput¬ 
ing the golden-rule rate and therefore that such methods 
will be necessary for the accurate simulation of electron 
transfer. 

The method is valid for any system with a simple 
potential-energy landscape, where either only one instan¬ 
ton exists or there are multiple instantons which are well 
separated. However, all instanton approaches will fail 
when the environment is fluxional, as is the case for 
electron-transfer in solution. Nonetheless, a semiclassi¬ 
cal instanton analysis can often be used to better under¬ 
stand the approximations involved in related methods. 
An example is the link between adiabatic instantons and 
RPTST,^® which has been shown to be exact in the ab¬ 
sence of recrossing. ®®’^^ A related semiclassical study of 
electron-transfer pathways^® has also helped improve at¬ 
tempts to describe nonadiabatic dynamics.In the fol- 
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lowing paper,^ we show how Wolynes’ quantum instanton 
method"^^ is related to our approach. 

When computing the semiclassical Green’s functions, 
we only considered imaginary-time trajectories and have 
therefore ignored real-time effects. For this reason the in¬ 
stanton approach should be considered to be a (quantum) 
transition-state theory. In the study of adiabatic reac¬ 
tion rates, it has been shown that ring-polymer molecu¬ 
lar dynamics (RPMD)®*^’*^ generalizes RPTST and hence 
the adiabatic instanton.^® Taking inspiration from this, 
it may be possible to combine such TST instanton meth¬ 
ods with nonadiabatic RPMD®^’®® to compute recrossing 
effects. 

The semiclassical Green’s function formalism is not 
limited to treating golden-rule rate problems. In forth¬ 
coming publications we will show how similar approaches 
can be applied to the Marcus inverted regime and to com¬ 
pute the rate constant in the adiabatic, large A limit. 
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Appendix A: Properties of determinants 


We shall make frequent use of the following expansions 
of the determinant of a block matrix: 

= |D| |A-CD-^Bl = |A| Id - BA'^Cl , (Al) 


A B 
C D 


which are valid as long as the inverses of D or A exist. 

When certain blocks in the determinant are equal, and 
A, B and 0 are square matrices of the same size, C and 
0 are vectors of the same length and D is a scalar, the 
following simplifications are found by rotating the basis 
set: 
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and also 


A B 
B A 


|A-B||A + B|. 
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